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Linear response theory for open (infinite) systems leads to an expression for the current response 
which contains surface terms in addition to the usual bulk Kubo term. We show that this surface 
term vanishes identically if the correct order of limits is maintained in the derivation: the system 
size, L, must be taken to infinity before the adiabatic turn-on rate 5 of the perturbation is taken 
to zero. In contrast to recent claims this shows that linear response theory for open systems is 
gauge-invariant without modification of the continuity equation. We show that a simpler derivation 
r-{ ' of the Landauer-Biittiker equations may be obtained consistently from the bulk Kubo term, noting 

that surface terms arising here are non-vanishing because they involve the opposite order of limits, 
, 5-^0, then L — * oo. 
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\Q • The connection between resistance and scattering properties of an open phase-coherent system, pioneered by Lan- 
dauer and Biittiker , has proved to be an invaluable tool in understanding the wealth of new transport phenomena 
that are observed in mesoscopic systems Given the S matrix elements for scattering from subband a in lead n to 
C ' subband a' in lead m, one can relate the currents and voltages in the leads by = J2n dmn Vn, where the conductance 
, coefficients are determined for to ^ n by 
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{9mm follows from current conservation) . The original derivation of this multiprobc formula by Biittiker Q for T = 
'"O is independent of the presence of a homogeneous magnetic field, because the cancellation between group velocity and 
Ch , density of states on which it relies still occurs at B ^ 0. The simplicity and generality of Eq. made it desirable to 
^ ' establish whether it followed from standard linear response theory (LRT) applied to open systems |Q. For the case 
of multilead structures in a non-zero magnetic field this was first shown by Baranger and Stone , and subsequently 
several different treatments were given 1^,0. 

Shortly thereafter, Sols |^ challenged the validity of these derivations, arguing that due to a boundary term arising 
only for open systems^ the Kubo formula could not be derived consistently without inserting an additional term into 
' the continuity equation. Indeed it was argued that this surface term violated gauge-invariance (since it vanished in 
some gauges and supposedly did not in others) and that the more general requirement of gauge-invariance necessitates 
the proposed modification of the continuity equation Q . Actually Appendix D of Ref. presented a proof that the 
surface term in question vanished identically in the non-interacting limit. The inconsistency according to Sols arose 
because subsequent manipulations leading to the Biittiker equations give rise to surface terms which were found to 
be non-zero and yet apparently had the same form as those shown to vanish in Appendix D. We have reexamined 
this issue and the closely related question of the invariance of LRT under global gauge transformations. Our results 
resolve the apparent contradiction by clarifying that the two limits which arise for LRT in open systems, the limit of 
the infinite system size and the adiabatic turn-on of the perturbation do not commute. The two apparently identical 
surface terms noted by Sols differ precisely in the interchange of these limits and wc show that they vanish in one 
order and not in the other as found in Ref. Hence Eq. (|l|) follows rigorously from LRT; a minor modification 
of our derivation proves the gauge-invariance of LRT without modification of the continuity equation. In addition 
our analysis leads to a derivation of Eq. (|l|) which is somewhat simpler than that in Ref. Q since it only employs 
properties of the scattering wavefunctions and makes clear the role of the adiabatic turn-on in breaking time-reversal 
symmetry. This derivation is similar to that given in Ref. Q except that we identify the central role of the above 
limits, clarifying in this way both the conceptual and technical steps leading to Eq. (nl). 

In Refs. an open conductor is represented by a finite sample region connected to N straight semi- infinite leads 

on which the external potential assumes constant but unequal values. We start with an adiabatically switched-on 
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DC perturbation of the form 



V{t) = e / p{r) 0(r) e** dr 



(2) 



where p(r) is the electron density operator. The expectation value of the current density caused by the perturbation 
in linear order can be written as 18(1 
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J(r) = -elim Jdt'e^^' Jd\ Jdr' {i{r, ~ihX) p{r' ,t')) 



(3) 



Applying the continuity equation ep{r, t') = — Vj(r, t') in the interaction picture and integrating by parts, one arrives 
at 
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3{r) = liin Jdt'e^'' JdX li{r,-ihX) - Jdr'i{r' ,t') ■ V (t){r') + J ds' j±{s' ,t') <j3{s') 



(4) 



where j±_ is the normal component of the current density, and the surface s' must be taken to infinity because the 
original integral, Eq. (|^), is of infinite range for the open system. However in the bulk term, which we denote by 
J^(r) because it leads to the Kubo formula, the spatial integral only has to be extended over the finite region A 
where the perturbing electric field E = — V'(/> is nonzero. To bring out the issue very clearly, let us now integrate by 
parts in the bulk term J^(r), noting that A is finite so the surface 6 A is not to be taken to infinity. We obtain 



J^(r) = lim Jdt'e^*' JdX /i{r,-ihX) -e Ja 



-e Jdr'p{r',t')cp{r',t') - Jds' j^{s' ,t') cp{s\t'\ 

A dA 



(5) 



The surface term in this equation differs from the one in Eq. ^ only in sign and in the absence of the s' ^ oo limit. 
Now we can determine the total current in lead m by integrating J(r) over a cross section of lead m. Let /^and 
/'„ denote the contributions of 3^ and the surface term in Eq. (0) to . The subtle point emerges once we take the 
limit 5 — > in the above expression for as the property V • — and the nonlocal Onsager relation for aoc in 
a magnetic field imply that the bulk term in Eq. (|5|) vanishes upon integration over a lead cross section, leaving only 
the surface term This is physically reasonable because the DC currents in the leads are uniquely determined by 
the voltages outside of A. But since the only restriction on A was that E = outside of it, nothing now prevents us 
from taking dA — >• oo. Thus becomes identical to — except for the interchange of the limits. If these two limits 
commuted we would obtain the absurd result Im = Im + = ^ • 

Note first that the bulk term in only vanishes when we take the limit (5 — > 0, so the order of limits leading to 
the surface expression for is uniquely determined. Conversely, the limit s' — > oo in Eq. must be taken before 
d can go to zero. The reason is that in calculating the LRT current as a trace over the unperturbed eigenstates, we 
assume that the system was in equilibrium in the distant past, which is achieved by (5 > 0. Matrix elements of the 
perturbation, Eq. (Q), in the unperturbed basis give rise to the spatial integral of infinite range. This is the asymptotic 
limit, which must consequently precede the adiabatic limit, S —>■ 0. It is not immediately obvious that this subtle 
difference really matters; therefore, we now evaluate both surface terms carefully, to show that vanishes while 
does not. Inserting many-body eigenstates \a> and executing the integrals over t' and A, one obtains from Eq. 



l!^ = J^(r) - ^h [dy„, ]™E ^ 

^ a,l3 
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Ea Ep- Ea+ ihS 



ds' </?|j(r)|a><a|jl(s')|/3> (/-(s'), 



(6) 
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where the integral extends over a cross section of lead m. The expression for /' is the same, only with dA replaced 
by s' — > oo. 

In the non-interacting limit, the statistical weights Pa can be replaced by Fermi functions f(Ea), and \a> become 
single-particle scattering states which, in lead Z, have the asymptotic form 



(7) 



Here, the label a consists of the energy E, subband index a and lead p of the incident wave. An analogous definition 
holds for -0^ = ipE'bq- "^^^ quantum wire eigenfunctions have the form 
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dk 
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(8) 



where x±a the transverse wavefunctions and kj_^ is the outgoing (+) or incoming (— ) wavenumber. With the 
above choice of normahzation, the symboUc sums over collective labels have the explicit form — » JdE. In 

the following, subscripts a imply a dependence on the energy variable E, while indices b belong to energies £". As 
was shown in Ref. current conservation implies the following properties for two wire eigenfunctions at the same 
energy in the asymptotic region: B] 



dyi (dj-L(r)IC'a') = S^'^Sa'a 



(9) 



with the notation — ±, and round brackets for the single particle current matrix elements. If ejj is the propagation 
threshold of subband a in lead I, Eq. (||) leads to integrals of the type 
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.f{E')-f{E) Fs{E',E) 



E' -E E' - E + m 



(10) 



where 



Fs{E',E) ^ |ds'(^^,,,lj^(r)|V^^,p)(V^^,^lj^(s')|V^'.,) 
s 



(11) 



contains the current matrix elements and thus vanishes at threshold. The surface S of integration is dA for I^, and 
is taken to infinity for /'. Although we wish to let dA oo eventually, we can do so only after (5 — > 0. Note also that 

f{E')-f{E) 



E'-E 



is non-singular on the real axis. 



We now evaluate the integral over E' in Eq. (10) by Cauchy's theorem. In general the integrand in Eq. (10) may 
have a complicated singularity structure in the complex plane, arising both from the "Matsubara" poles of the first 
factor and the singularities of the S-matrix contained in Fs{E, E'). However since asymptotically the wavefunctions 
have a plane- wave dependence on x' (the longitudinal coordinate in the lead), we shall see that contributions from 
the return contour may always be made to vanish as long as the correct half-plane is chosen to close the contour. 
This is so even if the return contour is a finite (not infinite) distance from the real axis |p^ . Since the locations 
of all the singularities of the integrand in Eq. (noh are independent of both the asymptotic and adiabatic limits 
except for the pole at E' = E — iTiS, for convenience we choose contours Ci, C2 (shown schematically in Fig. 1) 
which by assumption enclose only the latter pole. Hence Ci will contain no singularities at all, and terms which in 
the asymptotic limit may be closed in the upper-half plane will give zero contribution. C2 will enclose the pole at 
E' = E — iT%6 unless it approaches the real axis outside the range of the E' integration (i.e. unless E is less than the 
sub-band threshold relevant for the term in question) . By this convenient choice of contours we need only evaluate 
the residue at E' = E — ihS and see how it depends explicitly on the order of the asymptotic vs. adiabatic limits. 

The residue theorem yields with v = 1,2 
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f{E')~f{E) Fs{E',E) 



(12) 



E'-E E' -E + m 

-2.^ <5., Q{E-el) /(^ ' m - f{E) _ ^^^^ 

—ino 

where the step function Q enters since C2 encloses a pole only it E > e^, and (5^ 2 reflects the fact that no poles are 
enclosed by Ci. The additional minus sign is due to the negative sense in which C2 encloses the pole. Along C^, 
the limit 5 — > 00 enables us to apply the asymptotic forms of i/'a and ipfs in the second current matrix element in 
Fs{E',E) to get in lead I 



(V'Ub±(s')l^^'6,) = 5,i5j>i {e-a\3i-\tb)+5piY,Si,Mta\ji-\i+b') 

b' 



(13) 
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irrespective of whether (5 = or not in Eq. (|T^). For E' in the complex plane, the imaginary parts of the wavenumbers 
satisfy 

sgn(lmk^t^(i;')) = ±sgn(Im£"). (14) 



This can easily be checked for the special case kj_^ oc ±y/E' — e^, but holds even in asymmetric leads at i? 7^ where 
E' = el can occur for nonzero k. Note that the branch point at causes no problems since Fg{E',E) vanishes 

there. The first and third term of Eq. ( p^ depend on ^'^j, oc e*-''^ , so that they acquire an exponential decay 

factor e~'^^^ -'^ when integrated along the contour C2. For the other two terms in Eq. (|T^), we choose Ci to obtain 
positive imaginary parts in k'^;, (_£'). The asymptotic limit x' — * 00 causes the contributions from Ci and C2 to vanish 
identically. The result is that the real axis integrals over E' in Eq. ( |l0| ) vanish due to Cauchy's theorem unless C2 
encloses a pole. Its residue in Eq. (^2|) contains a factor e*-*'^^"*'''^'^ . For E ^ e^ we can expand the outgoing 
wavenumber as k'_^(£') — i5v^^{E) where the group velocity v{E) — j-^ is negative. We thus pick up a decaying 
exponential of the form 

FsiE-ihS,E) ocexp[-^ \v^\E)\ x'] . (15) 

This exhibits clearly the noncommutability of the limiting procedures: In /'„ we take a;' — > 00 first, causing Eq. ( p^ ) 
and with it I^^ itself to vanish. Moreover, since all boundary terms encountered in a gauge transformation Isl contain 
Eq. ( p^ ) with the above order of limits, this observation establishes the gauge invariance of LRT. 
On the other hand, I,^ is nonzero because now — > is performed first in Eq. (p^, yielding 

DO 



5^0 J E' -E E' - E + ih6 

4 

= -2m5,aQ{E - 4) f {E) Fs{E,E) (16) 

This resuh leads to Eq. (0), as we now show by taking the asymptotic limit S" — s- 00. The integral over s' in Eq. (|^) 
becomes a sum of integrals over isolated quantum wire cross sections, X^n/'^!/"' along which (p assumes its constant 
values Vn ■ Our expression for the total current at this stage already has the form Im = = 9mn Vn ■ The 
product of current matrix elements at the same energy in Fs{E, E) can be evaluated with Eq. ( |l3| ) if we let r in Eq. 
(^ go to infinity as well. Considering only n 7^ m, 7 of the 16 terms in Fs{E, E) vanish identically due to Kronecker 
delta factors like SpmSpn- When the sums over the incident wave parameters p, a or g, 6 are performed as prescribed 
by Eq. (^, five more terms drop out due to the unitarity of the S matrix, 

^ ^ ^ ^ ^np,a' c^mp.ac — Smn^a'a = 0. (17) 
p—1 c 

Here, we sum over all subbands propagating at energy E. This leaves the following four terms in Fs{E, E): 

5qm.5pmSnm,a'a^nm,b'b (C-bljilC-a) (?+a' IJ-L |C+b' ) + 
SqmSpnSmn,a'aSnm,b'b b'i IC+a' ) {C-a\j ±\C+b') + 

^pm5qnSnm,a'a^mn,b'b {^+a'\j-i-\C-b) {C+b' \j ±\^-a) + 

^qn^pn^rnn.b'b^^n.a' a {C+b'\j±\C'^a'){C\\j±\Cb), (18) 

where a sum over propagating subbands a', b' must be performed. The first two terms require integration along Ci 
in Eq. (|6|) and hence yield no residue. For the third term, Eq. (H) yields zero. In Ref. it is not noticed that this 
term does survive the contour integration and only vanishes due to Eq. (|^), which is not equivalent to the unitarity 
of the S matrix, but follows instead from the translational invariance of the asymptotic wire Hamiltonians. That this 
property of the model is important in deriving Eq. (^ was recognized solely in Ref. Q. In the present treatment, the 
remarks below Eq. (|lj) allow us to conclude that no other symmetries of the leads are essential. 

Our expression for now contains only the last term in Eq. (|l^). The conductance coefficients for m ^ n are 
then given by 

9mn=in dy„, dy,, (-27ri) dE Q{E - el) f (E) (19) 

•' •' p,q=la,b=l p 
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a'b' 

oo 

= jJdE [-f'{E)]Y,''\Sn.n,a'a\\ 

c-n o,'a 
^1 

where Eq. (^) has been used. This is the desired result, Eq. (|l|). Comparing the first and last term in Eq. (p^), we see 
that they differ only in the order of m and n. It is the sign of 5 that caused the first term to make no contribution 
to Eq. so that Qmn depends only on scattering probabilities from n to m, and not vice versa, as required by 

causality. Finally we comment that the vanishing of can also be demonstrated in the presence of interactions: the 
essential point is the smoothness of the integrand in Eq. (^0|) on the real axis for ^ ^ 0, which is preserved in the 
interacting case. 
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FIG. 1. The two possible contours that yield exponential 
decay in x' off the real axis for the various terms in Eq. (p^. 
Open dots: two possible locations of the pole at E' = E — ihS 
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